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Summary 


A new approach in estimating the Curie point Isotherm is 
developed by using the classical Gauss method o£ inverting a 
system of nonlinear equations. The method, slightly modified 
by the differential correction technique as described by Won 
(1981), will be used for a direct inversion of filtered Magsat 
data to calculate the crustal structure above the Curie depth, 
which is modelled as a magnetized layer of varying thickness 
and susceptibility. Since the depth below C:he layer is assumed 
to be non-magnetic, the bottom of the layer is Interpreted as the 
Curie depth. The method, once fully developed, tested, and com- 
pared with previous work by others, will be applied to a portion 
of the eastern U.S. when enough amount of Magsat data accumulate 
in the region. 


Review 

The basic objective in relating the aeromagnetlc field data with the 
structure of the Curie point isotherm is to compute the lower depth limit 
of magnetized masses in the earth's crust. Rocks lose their magnetism at 
the Curie temperature at which ferrlmagnetlc rocks become paramagnetic, and 
their ability to produce detectable magnetization dissappears. Thus, the 
deepest level in the crust containing materials with discernible magnetiza- 
tion is generally interpreted as the depth to the Curie point Isotherm. 

The Curie point is about 580“ C for magnetite. With appropriate titanium 
substitutions, Buddlngton and Ltndsley (1964) calculated an average Curie 
point ranging between 520“ C and 560“ C for rocks in the deep crust. It Is 
generally believed that the amount of geothermal heatflow should correlate 
with the Curie depth and thus, in turn, to the crustal magnetic field. 

Our main goal is, therefore, to determine the bottom shape of the 
magnetized crust from a magnetic anomaly map. Since the magnetic anomalies 
attributable to the bottom geometry are usually quite smaller and have much 
longer wavelengths than those produced by shallow geological variations, the 
problem is comparable to searching for a needle in a haystack. Early studies 
include those by Vacquier and Affleck (1941) and Bhattacharyya and Morley 
(1965). In both cases, each Isolated anomaly was filtered and separately in- 
terpreted by the empirical graphic method using a vertical-sided prism. 

A more sophisticated method was proposed by Bhattacharyya and Leu (1975s, 
1975b) . Their method requires an extensive initial filtering of the aero- 
magnetic data in both regional and short wavelength domains. The filtered 
data is subsequently divided into a large number of blocks. For each block, 
a two-dimensional spectrum and its moments are computed and compared with a 
model of an Isolated vertical-sided prism within a block in order to locate 
the corners of the body. The total amount of computation is tremendous since 
the method requires a two-dimensional Fourier transform for each block. Ap- 
plying the method to the Yellowstone National Park area, they produced the 
Curie isotherm map well correlated with the known geothermal area. 

Employing a similar technique, Shuey et al. (1977) concluded that it 
is essentially impossible to determine Curie depth with any resolution at 
all by fitting a vertical prism to a single anomaly. The Curie depths they 
derived could be changed by as much as 10 Km without violating the observed 
data. This conclusion is seemingly in conflict with those of Bhattacharyya 
and Leu (1975b) . 
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All methods reviewed here are commonly based on the assumption that 
there exl’^ts an Isolated magnetic source £or each anomaly. Each Individual 
anomaly Is assumed to be caused by a single vertical-sided prism (Bhattacharyya 
and Leu, 1975a, 1975b) or a truncated vertical cone (Shuey et.al. ,1977) . Such 
Isolated models are apt to generate spurious anomalies, particularly due to 
their unrealistically well-defined corners and vertical surfaces. These spur- 
ious anomalies can induce significant errors In either direct-modelling or 
spectrum calculation. 

Rock formations causing long wavelength magnetic anomalies at a depth 
close to the Curie point are more likely to have a continuous lateral distri- 
bution rather than Isolated blocks of well-defined geometrical bodies. A 
realistic model at this depth should manifest a continuous lateral^ distribution 
of magnetic materials having variable thicknesses and susceptibilities. 

Fluctuations in long wavelength magnetic anomalies can be attributed to 
lateral variations either In magnetization strength or in Curie depth. These 
double uncertainties make the task of simultaneously determining both the 
magnetizations and the Ci-srle depths very difficult, if not Impossible. Slmillar 
uncertainties apply to mrtsiy’ geophysical modelling theories, e.g., a thin magnetic 
dike for which the anomaly Is the same as long as the product of thlcknef/s and 
susceptibility remains the same. However, it can be shown that the statement Is 
no longer true If the dike has a considerable thickness for which case both the 
thickness and the susceptibility can be independently determined from observed 
data (Won, 1981). The following treatment is based on the classical Gauss 
method for solving non-linear equations (Corbato, 1965; Johnson, 1969) modified 
for digital total field magnetic data (Won, 1981) In order to Invert regional 
aeromagnetic data for a continuous crustal mass of varying thicknesses and 
susceptibilities . 


Proposed Method 

1. Mathematical Model for the Curie Depth 

Figure 1 shows the model which will be used for inverting aeromagnetic 
data. The model consists of laminated thick vertical prisms having flat top 
surfaces and linearly* connected Inclined bottom surfaces. The magnetic 
susceptibility below the lower boundary Is assumed to be zero so that the 
bottom geometry represents the Curie isotherm topography. Although data 
will be confined within the laminated block region, two semi-infinite slabs 
are added on either side in order to reduce the edge effects of the first 
and last blocks. The unknown parameters to be determined are the depth (h's) 
s.t each nodal point and the magnetic susceptibility (k's) of each prismatic 
body. 

The model is two-dimensional with an arbitrary strike angle with respect 
to the magnetic north. Data are assumed to be obtained at a constant alti- 
tude along a traverse perpendicular to the strike. Since the method uses 
total field aeromagnetic data, there Is no need for reducing the data to the 
polar anomalies . 

2. Mathematical Analysis 

Won (1981) formulated an inversion method for an Inclined dike of 
arbitrary thickness, depth of burial, dip angle, and susceptibility by 


employing Gauss’ method, as well as the differential correction technique. 
Starting from a set of initially guessed parameters, the inversion method 
self-corrects these parameters iteratively until the observed data fits 
best to the model. 

By combining two inclined dikes, we can derive a total field magnetic 
anomaly formula for a single vertical block having a flat top and an in- 
clined bottom. After some involving manipulation, we can derive the 
following formula for the entire laminated prismatic bodies as shown 
in Figure 1: 


F(X) 

2F, 


N-l 
S k, 

J-1 J 


[Aq^ + Bpj - s|{A-BAtj)q^ + (AAtj+B'< P|] 


+ k^ (a5q - Bpg) + (Aq^j + BPjj) 


( 1 ) 


F(X) denotes the total magnetic field anomaly measured at a distance X 
with respect to the origin. The first terra is the total contribution of 
N laminated bodies, while the second and third terms represent the con- 
tribution of the left semi-infinite slab and the right semi-infinite slab, 
respectively. Table 1 shows the mathematical notations used in equation 
(1). The formula can be used for a forward calculation of the total field 
anomaly fo,t any given set of variations in depths and susceptibilities. 

,3. Inversion of Anomaly by Gauss* Method and Differential Correction Technique. 

The unknown parameters include the depth profile (h's) and the lateral 
magnetic variation (k's) which are to be determined from a given profile 
obtained along a line perpendicular to the strike. For N laminated blocks, 
there are (2N-1) unknowns. All other quantities appearing in equation (1) 
are considered known. 


Suppose that the approximate values for h's and k's are assumed or 
otherwise established, even though the values may be significantly in error. 
Following a similar approach by Won (1981) , we shall describe a method of 
simultaneously improving all unknown parameters by using a technique com- 
bining the differential correction and the non-linear least-squares method. 


Let us write F(X) in equation (1) as 


F(X) *■ F (Xj h2» .... bjj} k£, .... (2) 

If we now chage hj (j * 1, 2, N) and kj (j * 1, 2, ..... N-l) by 

small amounts Ahj'^and Akj, respectively, the magnetic anomaly F(X) will 
change by AF such that 
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For M points of digitized total field magnetic data (1 
M >. 2N - 1), we attempt to minimize the quantity S 
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( 4 ) 


M 

S - Z (F. + AF. - G )' 
1-1 ^ 11 


with respect to Ahj and Akj. Differentiating with respect to Ah. and 
Akj» we obtain a oystem of (2N-1) simultaneous linear equations;^ 
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(5) 


The normal matrix in the first bracket has an order of (2N-1) and is 
symmetric and positive-definite. Once we solve equation (5) for the 
unknown differential correction vector CAh, , Ak, ], we obtain the improved 

values hT and such that ^ ^ 

J J ’ 


hj - hj + Ahj. 

kC - k, + Ah, 

J J 3' 


<6A) 

(6B) 


By replacing hj and kj with h| and k|, respectively, we can iterate the 
process until the RMS error is minimized. 

The expressions for derivatives in equation (5) are rather involving. 
Let us first define 




, ■ ly [<2AhjA + + (2AhjB - UjA)pj 

+ {(A - At^B)h^ - (AtjA+B)Xj}] 


(7A) 


and 


2 2 


[<2Ahj.iA + + (2Ahj.^B - Uj.jA)p'.i 

+ ^ ((A - - (AtjA+B)Xj) j . 

Derivatives with respect to h^ can now be written as 

aij “ “ *‘j-l j • 2, 3, (N-1), 

K 

“ j - 1. 

’ “Vl K for j - N . 


(7B) 
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Similarly, derlvatlvea with respect to kj*s can be written as 
^ 1^ . + Bpj - {(A - + (AAtj -f B)pj) 


S Tj ; for j - 2, 3, ,N - 1, (9A) 

- + (Aqg - Bp^); for j - X, (9B) 


The method described thus far does not necessarily converge to a correct 
final result; equation (3) is based on the assumption that the values of Ah. 's 
and Akj’s are small compared with the current values cf hj's and kj*s. There- 
fore when the current hj’s and k^’s are grossly wrong, equation (5) will produce 
large Ahj's and Akj's which may violate the assumption of differential correction 
and, consequently,^ the process may not converge. The safest way to avoid such 
a problem is to impose a certain allowable limit on each parameter. Such limit- 
ations may include; i) all depths are positive and less than, say, 100 Km, and 
il) all susceptibilities must be within, say, 0 and 0.01 c.g.s.. The conditions 
can even be tightened with more geological knowledge of the area. If any hj's 
or kj's go beyond its constraint, it may be allowed to change up to that limit. 

Unlike the Inversion approaches of individually correcting each unknown 
parameter by trial-and-error (e.g., McGrath and Hood, 1970) where the conver- 
gence to the absolute minimum error state depends on the cl\oice of Initial 
parameters, the present method simultaneously reduces all (2N-1) parameters 
regardless of their initial values. 

The method once developed will be first tested with the same data in the 
Yellowstone National Park Area used by Bhattacharyya and Leu (1975b) and the 
results will be compared critically. 


- 5 





Table 1 * 


Mathematical Notations and Abbreviations 


NOTATIONS; 


N 


d 

w 

F, 

I 

3 

X 


total number of nodal points 

Curie depth at j-th nodal point 

magnetic susceptibility of a block between J-th 
and (J+l)-th nodal points 

altitude of data (constant) 

vfidth of block (constant) 

the main Earth field In gamma 

Inclination of the main Earth field 

strike angle measured CCW from the magnetic north 

distance measured from origin 


ABBREVIATIONS ; 

A ■ sin 3 sin 21 

Xj - X - (j-l)w 

Uj ■ (w^-Ah"')/w 

s| ■ w^/(w^ + Ahj) 

rj - X? + d* 
j J 

io “ 

-p. - tan-l i - tan-1 \ 
*^0 x^ 


B ■ sin^ 3 cos^ I - sin^ I 


Ah, 


At, 


.^2 





h . + d 
J 

''j " ^j+1 

Ah./w 

j 

X ^ -f h"*^ 

’'j j 

ln(rj^l/rj) 

ln(r^/rN) 


tan 


-1 d 


*N 


tan 


“ 1 !!n 

*N 


tan”^ tan ^ •— — 

j+1 

_i - x.Ah, 

tan ‘aT'H? — H — ^ - tan . . 

Ahjhj + XjW h,h 


-1 


j''j+k ’‘j+l" 
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